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We-preseirt  an  algebraic  analysis  of  some  domain  decomposition  preconditioners  on  irregular 
regions.  We  analyze  a  preconditioner  proposed  in -[3]  for  the  interface  system  and  prove  that,  for  all 
L-shaped  regions  and  some  C-shaped  regions,  it  produces  a  convergence  rate  that  is  independent  of 
the  size  of  the  discretization  and  the  relative  shape  of  the  subdomains  (aspect  ratios).  Specifically, 
we  prov»  that  the  condition  number  of  the  preconditioned  capacitance  system  is  bounded  by  2.16 
for  all  L-shaped  domains.  We  ai-so  give  some  results  for  other  simple  irregular  geometries. 
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1.  Introduction 

We  consider  the  problem  of  solving  an  elliptic  partial  differential  equation  on  a  domain  that 
is  broken  up  into  rectangular  subregions.  By  using  domain  decomposition  or  substructuring  tech¬ 
niques,  the  problem  is  reduced  to  separately  solving  approximate  problems  in  the  subdomains  and 
updating  the  solution  at  the  interfaces  between  two  or  more  subregions.  For  the  class  of  domain 
decomposition  methods  considered  in  this  paper,  the  basic  idea  consists  of  the  following:  the  differ¬ 
ential  operator  is  discretized  on  a  grid  imposed  over  the  domain,  which  is  partitioned  into  several 
subregions.  Then,  by  applying  block  elimination  to  the  discretized  equations,  a  system  is  derived 
for  the  unknowns  on  the  interfaces  between  subregions.  This  system  is  sometimes  called  the  ca¬ 
pacitance  system.  Forming  the  right  hand  side  for  the  interface  system  requires  the  solution  of 
independent  elliptic  problems  on  the  subdomains.  For  certain  constant  coefficient  problems  on  reg¬ 
ular  domains,  fast  direct  methods  can  be  applied  to  the  solution  of  the  interface  system  [3,  4].  Such 
is  not  the  case,  however,  for  more  general  operators  or  irregular  domains.  For  efficiency  reasons  the 
system  must  then  be  solved  by  iterative  methods,  such  as  the  preconditioned  conjugate  gradient 
method.  Once  the  solution  is  known  on  the  interfaces,  one  more  elliptic  problem  must  be  solved 
on  each  subdomain  with  the  computed  values  as  boundary  conditions. 

In  [3],  an  eigenvalue  decomposition  in  terms  of  Fourier  modes  is  given  for  the  capacitance  matrix 
for  the  case  of  the  Poisson  equation  on  a  rectangle  divided  into  two  strips.  This  decomposition  is 
described  in  section  2.  In  this  paper,  we  are  interested  in  the  analysis  of  this  decomposition,  which 
we  will  cell  Me,  as  a  preconditioner  on  irregular  domains  and  in  particular,  we  want  to  study 
the  dependency  of  the  convergence  rate  on  the  gridsize  and  the  shape  of  the  domain.  Many  of 
the  preconditioners,  when  applied  to  an  L-shaped  region,  have  convergence  rates  that  are  bounded 
independently  of  the  gridsize.  The  bound,  however,  depends  on  the  relative  aspect  ratios  of  the 
subdomains.  For  example,  all  of  the  preconditioners,  except  for  Me,  are  known  to  deteriorate  when 
one  of  the  subdomains  becomes  narrow.  In  section  3,  we  show  that,  if  we  use  Me  as  preconditioner 
for  the  capacitance  matrix  on  any  L-shaped  region,  the  preconditioned  matrix  has  a  condition 
number  that  is  bounded  by  2.16,  independently  of  gridsize  and  aspect  ratios.  Given  an  L-shaped 
region,  there  are  two  ways  of  separating  it  into  two  rectangular  subregions.  We  prove,  also  in 
section  3,  an  interesting  property  of  the  preconditioner  Me,  namely  that  the  convergence  rate  is 
not  affected  by  the  way  we  choose  to  subdivide  the  domain.  In  section  4,  we  discuss  the  extention 
of  some  of  the  results  in  section  3  to  C-shapes.  In  the  proofs  of  sections  3  and  4,  we  often  use  a 
common  operator,  which  describes  the  interaction  between  two  perpendicular  interior  interfaces. 
This  operator  is  analyzed  in  detail  in  the  appendix. 

2.  The  interface  operator  and  its  preconditioners 

In  order  to  illustrate  the  method,  we  will  apply  the  process  described  above  to  a  simple  region 
0.  which  can  be  decomposed  into  two  rectangles  fli  and  O2,  with  interface  T3,  as  shown  in  fig.  1 . 
Let  the  linear  system 

Au  =  f  (2.1) 

represent  the  discretization  of  the  differential  operator  on  O.  By  ordering  the  variables  in  f 1]  and 
Q2  first  and  then  those  in  r.3,  the  system  (2.1)  can  be  written  in  block  form  as: 

f  -111  ^13  \  /  «1  \  /  fl  \ 

[  A22  .423  )  u2  =  1/2  I  (2.2) 

\  -4  kj  .4 [3  .433/  w  W 


where  the  indexes  for  u  and  /  corespond  to  gridpoints  in  0-,,  0^  and  f  j,  respectively.  Based  on 


Figure  1:  The  domain  ft  and  its  partition 


the  following  block  decomposition  of  the  matrix  in  (2.2): 

A  = 


aJ3  c 

where  C  is  the  Schur  complement  of  T33  in  .4,  i.e. 


22^23  I 


4-1 


C  —  4.33  —  /I13  —  j423^22J  ^23 


T  A- 1 


the  system  (2.2)  can  be  solved  as  follows: 
Step  1:  Compute 


g  -  fa-  /1  —  A23A22  fa 


T  ,-l 


and  solve 


Step  2:  Solve 


and 


Cu3  =  g 


A\\U\  =  fi  -  A13U3 


(2.3) 

(2.4) 

(2.5) 

(2.6) 

(2.7) 


A22U2  —  f2~  A23U3  (2-8) 

The  computation  of  g  by  (2.5)  and  uj  and  «2  by  (2.7)  and  (2.8),  require  the  solution  of  independent 
problems  on  the  subdomains.  The  matrix  C  given  by  (2.4),  also  called  the  capacitance  matrix,  is 
dense  and  expensive  to  compute.  It  is  possible,  however,  to  compute  the  action  of  C  on  a  vector 
v  at  the  cost  of  solving  problems  on  the  subdomains  with  boundary  conditions  on  T  given  by  v. 
Therefore,  the  interface  system  (2.6)  is  often  solved  by  preconditioned  conjugate  gradients  (PCG). 
Since  each  iteration  involves  solving  problems  on  the  subdomains,  it  is  essential  to  keep  the  number 
of  iterations  low.  For  this  reason,  much  effort  has  been  devoted  recently  to  the  construction  of  good 
preconditioners  for  the  capacitance  matrix  [6,  1,  7,  3,  4].  Many  of  the  preconditioners  proposed  are 
spectrally  equivalent  to  the  exact  boundary  operator.  They  therefore  yield  convergence  rates  that 
are  bounded  independently  of  the  gridsize.  The  method  is  particularly  suited  to  problems  for  which 
the  subproblems  can  be  solved  efficiently,  for  example,  when  the  operator  has  separable  coefficients. 
When  the  subdomain  problems  cannot  be  solved  efficiently  but  they  can  be  approximated  by 
separable  operators,  it  is  possible  to  derive  block  preconditioners  for  the  original  system  based  on 
preconditioned  for  the  interface  system  [8,  2,  5]. 


2 


t 


i 


In  [3j,  the  case  of  a  constant  coefficient  operator  on  a  rectangular  domain  divided  into  two  strips 
is  analyzed.  For  this  simple  case,  it  is  shown  that,  for  many  of  the  preconditioners  proposed  in  the 
literature,  while  the  condition  number  of  the  preconditioned  system  can  be  bounded  independently 
of  the  gridsize  h  for  a  fixed  domain,  it  can  grow  as  a  function  of  the  aspect  ratio  of  the  subdomains. 
Roughly  speaking,  the  aspect  ratio  of  a  rectangle  is  the  ratio  between  its  height  and  its  width  (note: 
for  one  of  the  preconditioners  proposed  in  [l],  the  bound  grows  when  only  one  of  the  subdomains 
becomes  narrow).  A  fast  direct  solver  for  C  based  on  Fourier  analysis  can  be  derived  from  the 
exact  eigenvalue  decomposition  of  the  capacitance  matrix.  This  operator  takes  aspect  ratios  into 
account  and  solves  exactly  the  interface  problem  for  the  case  of  constant  coefficients  on  a  rectangle 
divided  into  t%vo  strips.  It  is  therefore  proposed  in  [3]  to  apply  it  as  a  preconditioner  for  interface 
systems  on  irregular  regions  or  for  variable  coefficient  operators.  We  will  call  this  preconditioner 
Mr.  For  the  case  of  a  five  point  finite  differences  discretization  of  the  Poisson  equation: 

~UXX  —  Uyy  —  f  (2-9) 

on  a  regular  grid  of  size  h  =  Me  has  a  decomposition  of  the  form: 

Mc  =  1V„  A  Wl  ,  (2.10) 

where  A  is  a  diagonal  matrix  and  Wn  is  the  matrix  of  sine  modes  of  dimension  n,  whose  elements 

are  given  by:  _ 

I  2  .  ijir 

Wij  =  V  nTI  Sm  nTT  (2-U) 

for  i,  j  =  1, . . . ,  n. 

Given  integers  n,  mi  and  m2,  define 


(i  +  iTl+l  1  +  7?2+l \  I 

A,(n,  mi,  m2)  =  — pr  +  !  1  7m2-+f  j  ^  +  X  (2’12) 

where 

<’'=W(>TT]f)  (2'13) 

and 

v  =  (*  +  f -\[°i  +  r)  ■  <214> 

The  eigenvalues  of  Me  are  given  by 

A  j  —  A_,(n,  mi,  rri2) 

for  j  =  1, . . .,  n,  where  mi  and  m2  are  the  number  of  grid  points  in  the  {/-direction  in  fti  and 
respectively. 

The  preconditioners  proposed  in  [6]  and  [7]  have  the  same  eigenvectors  as  (2.10),  but  the 
eigenvalues  are  those  of  the  square  root  of  the  one-dimensional  discrete  Laplacian,  namely  v/cr“  in 

[0]  and  \f <7  j  +  in  [7],  For  the  case  of  the  Poisson  equation  (2.9),  it  can  be  proved  that  one  of 
the  preconditioners  given  in  [1]  also  Iras  a  decomposition  of  the  form  (2.10).  The  eigenvalues  A;  for 
this  operator  can  be  obtained  by  setting  m2  -  m.i  in  (2.12),  i.e.  Aj(n.  m,,m  1 ).  This  preconditioner 
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is  therefore  exact  for  the  case  of  a  rectangle  divided  symmetrically  into  two  identical  rectangular 
subdomains. 

3.  L-shaped  regions 

In  this  section,  we  describe  the  interface  operator  and  its  preconditioners  for  an  L-shaped 
domain,  the  simplest  irregular  shape  that  can  be  decomposed  in  rectangular  subregions.  Consider 
the  Poisson  equation  on  the  region  fi  of  fig. 2.  It  is  clear  that  either  interface,  1%  or  IV  will  divide  the 
domain  into  two  rectangles.  We  might  ask  ourselves  two  questions:  is  a  particular  decomposition 
better  than  the  other?  And  how  does  the  convergence  rate  depend  on  the  mesh  size  and  the  aspect 
ratios  of  the  subdomains?  We  will  show  that  for  the  particular  preconditioner  Me  we  analyze, 
the  two  decompositions  produce  iteration  matrices  with  the  same  convergence  rate.  We  also  give 
a  bound  for  the  condition  number  that  is  independent  of  the  mesh  size  and  the  subdomain  aspect 
ratios. 
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r4 
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fl2 

r5  n3 

Figure  2:  L-shaped  domain 

Let  the  linear  system 

Au  =  /  (3.1) 

represent  a  standard  second  order  five  point  discretization  of  the  differential  equation  on  a  regular 
grid  imposed  on  the  domain  f2.  Let  us  first  consider  the  domain  as  the  union  of  two  rectangles 
divided  by  the  interface  IV  An  interface  system  of  the  form 


C4u4  =  gA  (3.2) 

can  be  derived  for  the  variables  on  by  the  process  of  block  elimination,  similarly  to  equations 
(2.5)  to  (2.8). 

Similarly,  we  can  consider  the  domain  fi  as  the  union  of  two  rectangles  divided  by  the  interface 
rs  and  an  interface  system  of  the  form 


C5U5  =  g  5 


(3.3) 


can  be  derived  for  the  variables  on  IV 

On  the  other  hand,  by  reordering  the  gridpoints  on  the  subdomains  first  and  then  those  on 
the  interfaces  I\j  and  IV  A  can  be  written  in  block  form  as: 


A  = 


( 3.4) 


4 


where 


•In 


An  = 


M2 


-4  33 , 


-4p  =  (  A*4 


*55 


and  P  = 


,)('  -ir)  ■ 

(3.5) 

L,  i.e., 

M4  —  ^24^22  A25  A 

"A  25 .4  22*  '424  Ms  J 

(3.6) 

M4  —  A  24  A*>2  A24 

(3.7) 

A22  A25  —  aJ5  A33x  A35. 

(3.8) 

The  matrix  A  of  (3.4)  can  be  decomposed  as  follows: 

A=(pf  C45 

where  C45  is  the  Schur  complement  of  /Ip  in  A,  i.e., 

C45s  Ar-PTA~1P= 

with 

M4  —  .4  4  4  — 

and 

Ms  =  A55  -  A\sA 

The  matrix  M4  would  be  the  capacitance  matrix  for  T4  if  the  domain  M  were  absent.  Similarly, 
Ms  would  be  the  capacitance  matrix  for  Ts  if  the  domain  were  absent.  In  fact,  they  are  nothing 
but  the  preconditioner  Me  described  in  the  previous  section.  Both  A/4  and  M5  have  eigenvalue 
decompositions  of  the  form  (2.10).  According  to  the  definition  (2.12),  the  eigenvalues  of  M4  are 
given  by  A j{n.mx,m2)  for  j  =  l,...,n  and  its  eigenvectors,  by  Wn.  The  eigenvalues  of  Ms  are 
A, (m2,  n,  nz)  for  i  -  1, . . . ,  m2,  with  eigenvectors  given  by  Wm2. 

The  matrix  C4  of  (3.2)  is  the  Schur  complement  of  ,444  in  A,  but  it  can  also  be  written  as 
the  Schur  complement  of  M4  in  C45.  Similarly,  C5  is  the  Schur  complement  of  A55  in  .4,  but  it 
can  also  be  written  as  the  Schur  complement  of  M5  in  C45.  Therefore,  we  can  derive  the  following 
expressions  for  C4  and  Cs'- 

Lemma  3.1.  The  interface  matrix  for  T4  in  ft  can  be  written  as 

C4  =  M4  -  BT  Mi  lB  , 

where  B  =  A25A22  A?4.  Similarly,  the  interface  matrix  for  r5  in  ft  can  be  written  as 


(3.9) 


Cs  =  Ms  -  BM71  Bt 


(3.10) 


The  preconditioner  proposed  in  [3]  for  C4  in  (3.2)  would  correspond  to  Me  =  M4  and  similarly, 
Me  =  Ms  for  Cs  in  (3.3).  Since  Me  is  positive  definite,  we  can  choose  \fMc  as  a  symmetric 
preconditioner.  Let  us  define  the  preconditioned  matrices: 


By  (3.9).  we  have 


c4  =  m41/2c4m4  1/2 

and 

Cs  =  M-1/2CsMTl/2  . 

(3.11) 

C4  =  In  -  BTB 

and 

Cs  =  fm2  -  BBt  . 

(3.12) 

■j 


where 


b  =  m;'/2aZsa£a2*m;1/2  .  (3.13) 

If  we  choose  r4  as  the  interface,  at  each  iteration  subdomain  problems  will  be  solved  on  and 
fl-2  U  r5  U  $^3.  Similarly,  if  we  choose  Ts  as  the  interface,  at  each  iteration  subdomain  problems 
will  be  solved  on  U  r4  U  and  Clj.  The  work  per  iteration  is  therefore  comparable  for  both 
ways  of  splitting  the  domain.  We  will  next  show  that,  by  solving  (3.2)  with  preconditioner  .\/4 
and  (3.3)  with  preconditioner  A/5,  both  systems  are  also  equivalent  from  the  convergence  point  of 
view.  Therefore,  in  a  general  case,  there  is  no  a  priori  reason  to  prefer  one  way  of  decomposing  the 
domain  over  the  other. 

If  n  =  m2,  this  fact  is  not  surprising,  considering  that  both  interface  systems  have  the  same 
order  and  it  is  easy  to  see  that  C4  =  C5.  It  is  not  obvious,  however,  whether  one  way  of  decomposing 
the  domain  should  be  prefered  when  n  ^  m2.  As  it  turns  out,  even  in  this  case  the  asymptotic 
convergence  rate  is  the  same  for  both  systems,  because  the  matrix  C45  of  (3.6)  satisfies  the  following 
theorem: 

Theorem  3.1.  Consider  the  following  symmetric  positive  definite  (SPD)  system,  written  in  block 
form : 

d  ;)(:)■(!)  • 

where  the  blocks  ,4  and  B  are  square  matrices.  Also,  define  the  Schur  complement  systems: 

(A- BC~1BT)x  =  f  -  BC-'g  (3.15) 

and 

(C  -  BrA~1B)y  =  g- BTA~1f  .  (3.16) 

Consider  the  solution  of  (3.15)  by  the  following  fixed-point  iteration,  with  splitting  matrix  given 
by  ,4:  given  an  initial  guess  x° ,  define  the  i-th  iterate  as  the  solution  to: 

Ax'  —  f  —  BC~l g  +  BC~l  BT x'~J  (3.17) 


for  i  =  1,2,... 

Similarly,  given  y°,  define  the  i-th  iterate  of  a  fixed-point  iteration  for  solving  (3.16)  with 
splitting  matrix  given  by  C,  as: 

Cy'  =  g  -  BTA~lf  +  BTA~1By'-1  (3.13) 


for  1  =  1,2,... 

Then,  the  two  iterations  are  convergent.  Moreover,  they  are  equivalent  in  the  sense  that  for 
any  given  initial  guess  x°  for  (3.17),  there  exists  an  initial  guess  y°  for  (3. IS),  such  that  for  all 
i  =  0,  1, ...  we  have: 

y'  =  qy  +  Pyx'  ,  (3.19) 

where  qy  =  C~lg  and  Py  —  —C~1Br  and 

l|e'+1lU<IKHc  <||^|U  .  (3.20) 

where  e‘x  =  x‘  —  x  ,  ey  =  yl  —  y  and  [|  uj| ^4  denotes  the  A-norm  of  a  vector  u,  i.e.  VuT  A  u. 

Completely  analogous  results  also  hold  for  x'  and  elx,  given  an  initial  guess  y°  for  (3.16). 

Proof.  Given  x°,  define  y°  —  qy  +  Pyx° .  By  induction,  we  can  see  that  (3.19)  satisfies  (3.18)  for 
every  i  >  1. 
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From  the  classical  matrix  iterative  analysis  for  the  convergence  of  block  Gauss- Seidel  iteration 
for  SPD  matrices,  it  can  be  shown  that  the  two  iterations  converge.  Also,  since  the  matrix  of  (3.14) 
is  SPD.  so  are  the  blocks  .4  and  C  and  their  corresponding  Schur  complements.  Therefore,  A-1/2 
and  C-1/2  are  well  defined  and 


UA-i/2BC-i/2||2  <  x 

(3.21) 

We  can  also  prove  that 

Ae'*1  =  BC~l  BTe'x 

and 

J  —  p  J 
cy  —  1  ycx 

Therefore,  we  have 

.41/2ej.+1  =  -(A~1/2BC-l/2)Cl/2ei 

(3.22) 

and 

Cl!2e[  =  -(C~1/2  BT  A~1/2)A1,2eix  . 

(3.23) 

Using  (3.21),  (3.22)  and  (3.23),  we  can  prove  (3.20). 

1 

When  an  iterative  method  such  as  PCG  is  used,  the  rate  of  convergence  depends  on  the 
condition  number  of  the  corresponding  preconditioned  matrix  in  (3.11).  By  applying  the  last 
theorem  to  C45  and  by  using  (3.12),  we  can  conclude  the  following: 

Theorem  3.2.  Solving  both  systems  (3.2)  and  (3.3)  with  preconditioners  of  the  form  Me  produce 
equivalent  asymptotic  convergence  rates.  Moreover ,  by  (3.12),  we  have 

/C(C4)  <  - - 

'  ~  1  -  P^||2 

(3.24) 

and 

K(Ci)  <  --4-t  -  . 

l-PrS||2 

(3.25) 

Numerical  computations  show  that  the  singular  values  /?,-  of  B  decrease  very  quickly  with  the 
index  i.  Therefore,  in  practice,  only  a  few  eigenvalues  of  C4  and  C5  are  different  from  1,  which  leads 
to  rapid  convergence  of  the  PCG  method  when  applied  to  either  matrix.  For  example,  for  the  L- 
shaped  region  with  corners  at:  (0,  0),  (3, 0),  (3,  0.25),  ( 1, 0.25),  (1, 1.25)  and  (0, 1.25),  for  n  =  31  and 
63,  table  1  shows  the  singular  values  of  B  and  the  eigenvalues  of  C5,  computed  in  single  precision. 

Our  conclusion  is  that  either  way  of  decomposing  an  L-shaped  region  into  two  rectangles 
produces  the  same  convergence  rate,  when  preconditioner  Me  is  used.  Moreover,  we  will  be  able  to 
give  an  analytical  bound  on  the  condition  number  of  the  preconditioned  capacitance  matrix.  This 
bound  is  derived  from  a  bound  on  the  norm  of  the  operator  BT B. 

But  first,  we  will  give  an  expression  for  the  elements  of  a  unitary  transformation  of  B.  Let 
the  elements  of  the  matrix  Wn  be  given  by  (2.11)  and  similarly,  define  the  elements  of  fUm2  by 
replacing  n  by  m2  in  (2.11). 

The  operator 

Q\E  =  A  25^22' -^2 ‘I  - 

which  is  part  of  the  definition  (3.13)  of  B,  is  the  operator  that  takes  boundary  values  on  the 
interface  r.t.  solving  a  Poisson  problem  on  fi2  and  then  takes  the  values  of  the  solution  at  the 
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! 

!■ 

to 


n  =  3 1 ,  m2  =  7 

n  -  63,  m2  =  15 

sv  of  D  (r(Cs) 

sv  of  B  <t(  C5 ) 

0.18204  0.96686 
0.03868  0.99850 
0.00514  0.99997 
0.00045  0.99999 
0.00002  1.09000 
0.00000  1.00000 
0.00000  1.00000 

2.165E-01  0.95312 
6.816E-02  0.99535 
1.578E-02  0.99975 
2.971E-03  0.99999 
4.607E  04  0.99999 
5.863E-05  1.00000 
6.082E-06  1.00000 
5.093E-07  1.00000 
3.610E-08  1.00000 

Table  1:  Eigenvalues  of  preconditioned  capacitance  system 
for  an  L-shaped  region 

gridpoints  which  are  adjacent  to  r5.  It  is  possible  to  derive  the  elements  of  Qs'e  when  it  is  pre 
and  post-multiplied,  respectively,  by  the  matrices  Wm2  and  VFn.  The  elements  of 

Wm2QNEWn 


are  given  by 

-  2  sin  sin  n+T 

q'3  y(m2+l)(n+l)  a\n)  +  a,(m2> 


(3.26) 


for  i  -  1, . . . ,  m2  and  j  -  1, . . n.  A  proof  of  (3.26)  can  be  found  in  the  appendix  (see  lemma  5.2). 

For  any  given  integers  n,  mi  and  m2,  let  A(n,mi,m2)  be  defined  by  (2.12),  where  7;  is  given 
by  equation  (2.14).  By  using  (3.26),  it  is  easy  to  prove  the  following  lemma: 


Lemma  3.2,  Let 

V  =  WmjBWn  .  (3.27) 

Then ,  \\V\\2  =  |(5||2  and  the  elements  of  the  matrix  V  are  given  by 


sin 


7712+1 


sin 


JTX7T 

n+1 


y/{n+  l)(m2+  1) 


(<rjn)  +  Crt(m2)) 


(3.28) 


for  i  =  1, . . .,  m2  and  j  =  1, _ n,  where  s*4)  =  ^/[ Ay ( re,  mi,  m2)|  and  s-5)  =  v/|A,(m2,  n,  n3) 


As  equations  (3.24)  and  (3.25)  suggest,  in  order  to  find  a  bound  for  the  condition  number  of 
the  preconditioned  capacitance  system,  we  need  to  bound  the  norm  of  B,  or  V,  Since  we  have  an 
expression  for  the  elements  of  V,  we  can  bound  ||Vj|i  and  IjVIjco  and  then  use  the  property: 


v'ii2  <  vwunz 


The  results  are  summarized  in  the  next  theorem.  A  proof  can  be  found  in  appendix  B: 
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theorem  3.3.  Define  the  aspect  ratio  for  the  domain  fi2  in  fig.  2asa='^l.  Then. 

a)  ;|l'l|i  <  v/a  0.733  and  P'|U  <  0.733. 

b)  ||5Tfl||:  <  \\D\\]  =  IH'Hj  <  lll'IMIl'Hoo  <  0.54. 

c)  For  all  gridsizes  and  all  L-shaped  regions. 

AJ(C.|)  <  2.16  and  X(C5)  <  2.16  .  (3.29) 

I 

In  our  experiments  on  L-shaped  domains  with  many  different  aspect  ratios,  condition  numbers 
larger  than  1.2  have  not  been  observed.  The  bound  0.54  in  b),  however,  is  fairly  tight  for  |[V'||i  j|l '||x.. 
as  was  shown  by  numerical  experiments  with  large  values  of  n  and  m 2 .  Therefore,  if  a  tighter  bound 
is  desired  for  the  condition  number,  one  would  need  to  bound  the  2-norm  of  Dr B  directly. 

We  would  also  like  to  discuss  briefly  how  the  parameter  n 3  (or,  respectively,  m , )  affects  the 
performance  of  preconditioner  .V/t  (Ms).  Clearly,  as  713  tends  to  zero  for  large  m2,  the  domain  S3 
approaches  the  shape  of  a  perfect  rectangle.  The  preconditioner  M 4  should  reflect  this  by  becoming 
the  exact  boundary  operator.  In  other  words,  AC(C4)  should  approach  one.  We  can  verify  that  this 
is  the  case  as  follows:  r,_,  in  (3.28)  depends  on  773  only  through  A, (77? 2-n-na)  (defined  in  (2.12)). 
When  the  aspect  ratio  tends  to  zero  (i.e.  Q3  becomes  thinner),  A,(m2, 71,713)  tends  to  infinity 

and  therefore  utJ  tends  to  zero.  However,  we  can  see  that  this  dependency  is  very  weak,  because 
Aj(  77t 2 . 7?,  773)  tends  rapidly  to  an  asymptotic  value  independent  of  773  when  such  aspect  ratio  grows. 
Only  the  fact  that 

AJ(7n2,  n,  713)  >  2^/trJ  (3.30) 

is  used  in  the  proof  of  theorem  3.3,  which  is  true  for  all  n3.  The  discussion  above  implies  that  the 
performance  of  ,U4  as  a  preconditioner  for  C4  is  fairly  independent  on  how  irregular  the  region  is. 

Incidentally,  for  the  other  preconditioners  mentioned  in  this  paper  [6,  1,  7],  the  preconditioned 
capacitance  matrix  always  has  the  form  X  +  BT B,  for  some  operator  B  to  which  the  bounds  (a)  and 
(b)  of  theorem  3.3  can  also  be  applied,  as  long  as  (3.30)  holds.  The  bound  given  in  (c),  however, 
does  not  hold  for  other  preconditioners,  for  which  the  norm  of  X  may  grow  when  the  aspect  ratio 
a  of  the  domain  fl2  decreases  (see  [3]  for  an  example  on  a  T-shaped  region). 


4.  C-shaped  regions 

Some  of  the  expressions  and  results  of  the  previous  section  are  more  general  than  they  appear 
and  they  can  be  used  as  basic  components  for  more  complicated  regions  that  are  unions  of  rectangles. 
For  example,  a  C-shaped  region  can  be  subdivided  as  indicated  in  fig. 3. 

Similar  to  L-shaped  domains,  the  region  of  fig. 3  can  be  separated  in  three  rectangles  by  either 
r6  and  Tr,  or  Tg  and  Tg.  By  ordering  the  variables  in  ft;,i  <  5  first  and  then  those  on  Tj,  6  <  j  <  9. 
the  matrix  A  that  represents  the  discrete  differential  operator  on  can  be  written  in  block  form 
as  in  (3.4),  where 


•4,-  = 


-4  16  ‘4 18 

.4  26  .4  27 
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Figure  3:  C-shaped  domain 


A  system 


(4.2) 


can  be  derived  by  block  elimination  for  the  interfaces  and  T 7,  where  Cgr  is  the  Schur  complement 
in  .1  of  the  blocks  .lee  and  A-t-  In  [4],  a  multistrip  operator  M67  is  described,  which  solves,  exactly, 
the  problem  on  a  rectangle  divided  into  three  strips  (fi!,fi2  and  ft3).  We  will  analyze  Me-  as  a 
preconditioner  for  C&7.  The  operator  A/67  has  the  following  block  structure: 


A/67 


(4.3) 


where 


and 


#6  =  ^66  -  -4 16^7/  ^16  -  ^26‘^22  ^26  , 
—  T77  —  Aj7A3^  A37  —  AlrA.i  -427 

S  =  —  A^6A22  -427  • 


The  blocks  He,  H 7  and  5  have  eigenvalue  decompositions  of  the  form  (2.10),  with  eigenvalues  given 
by  \j(n,mi,  m2),  A j(n.  m 2,  m3)  and 


6}(n,  m2) 


(4.4) 


respectively,  for  j  =  1 . n.  (See  lemma  5.1  in  appendix  A). 

Similarly,  a  system 


can  be  derived  for  the  interfaces  Fg  and  Tg,  where  Cgg  is  the  Schur  complement  in  .4  of  the  blocks 
,4><g  and  Agg.  The  system  (4.5)  can  be  preconditioned  by  a  block  diagonal  preconditioner  A/go.  with 
diagonal  blocks  A/g  and  Mg.  .!/«  is  the  exact  interface  system  for  Tg  with  respect  to  the  subdomains 


1* 

f 


Sl[  and  fi.t,  and  .\/9  is  the  exact  interface  system  for  Ty  with  respect  to  the  subdomains  i >;j  and  ft5. 
Both  Ms  and  ,\/.j  have  decompositions  of  the  form  (2.10). 

It  can  be  easily  shown  that  Cp.  the  Schur  complement  of  the  blocks  .4p  in  .1.  ran  be  written 
in  block  form  as: 


C'r  - 


Ms  7 

Qse 

0 

0 

Q.we 

QIse  0 

Ms 

0 

0  QIe 

0 

Mg 

where  Qse  —  .dfg.lJ'/Ajg  and  Q,\e  =  Again,  by  applying  theorem  3.1.  we  have  that 

both  ways  of  dividing  the  domain  are  equivalent,  in  the  sense  that  initial  residuals  can  be  found  such 
that  the  same  number  of  iterations  are  necessary  when  PCG  is  applied  to  C&~  with  preconditioner 
.l/er  than  when  PCG  is  applied  to  Cg9  with  preconditioner  A/g9.  The  preconditioned  interface 
operator  for  Fg  and  r7, 

1/2 


C67  =  Mq71/2C67Mq-1/ 


can  be  written  in  the  form 
where  B  £  ft(mi+m3)x2n  an(j 

B  = 


-  t) 

C67  =  I  -  btb 


M&  0  Yl/2(QTSE 
0  Mg)  V  0 


QT 


0 

T 

NE 


A/e- 


1/2 


Similarly,  the  preconditioned  interface  operator  for  for  Tg  and  r9, 


C89 

can  be  written  in  the  form 


_  (Ms 


0  Mg 


0  ^ 

u9) 


-1/2 


Cs9 


Ms  0  \  ~1/2 

0  Mg 


Cgg  =  /  -  BB 7 


The  condition  numbers  of  C&7  and  Cg9  are  boundeu  by 

1 


and 


IC(C67)  < 


A1(Cs9)  < 


1  - 

1 


i-\\btb\\2 

Define  V  as  the  following  unitary  tranformation  of  B: 


V 


-( 


W, 


m  1 


0 


B 


\Vn 


0 


Then  lll'H 


0  wm. )  \  0  wn 

|Z?||.  The  matrix  V  can  be  written  as  a  block  two  bv  two  matrix 

1 r  _  I  I  66  ^  67 

~  Ire  Irr 


(4.6) 


(4.1 


(4.8) 


(4.9) 


4.10) 


1 1 


Figure  4:  Interaction  between  interfaces 

whose  block  elements  have  expressions  similar  to  the  matrix  V  for  L-shaped  regions,  namely, 

Vee  =  \VmiM;l,2QTSEWnR6 


Vsr  -  WmiM-1/2QTSEWnR _ 
Vre  =  WmiM~xl2QTNEWnR_ 
V77  =  WmzM~1/2QlEWnRr 


(4.11) 


where  R$,  R7  and  are  diagonal  matrices  such  that: 


R6  R- 
R-  R7 


:HWo 


Wn  0 
0  W„ 


For  the  case  when  mx  —  m3  <  m2,  a  simple  expression  can  be  found  for  R6,  R7  and  R . ,  namely 
R$  =  R7  =  R+,  with  the  diagonal  elements  of  R±  given  by 


rf  =  - 

J  9 


/*JTW 


where  Aj  is  A }(n,  mx,  m2),  given  by  (2.12)  and  Sj  is  6j(n,  m2),  given  by  (4.4).  Arguments  simi  ir  to 
those  in  theorem  3.3  can  be  applied  to  give  the  following: 

Theorem  4.1.  Consider  a  C-shaped  region  like  fig. 3,  where  mx  =  m3  <  m2  and  a  is  the  asx  ret 
ratio  for  the  domain  ftj  or  173  in  the  picture,  i.e.  a  =  2jL±!_  Then, 

&)  ||F||,  <  y/a  0.7877  and  UVH*,  <  ^5  0.7877  . 

b)  ||V^V||2  <  HVII*  =  ||£||*  <  H5IMI5IU  <  0.62. 

c)  IC(C6t)  <  2.63  and  /C(C89)  <  2.63  for  all  gridsizes  and  all  C-shaped  regions  such  that 
m  1  =  m3  <  m2. 

Proof.  In  appendix  B. 


5.  Appendix  A:  The  interaction  between  interior  edges 

In  this  appendix,  we  define  the  operators  that  represent  the  interaction  between  two  interfaces 
of  a  given  subdomain.  Consider  the  rectangular  region  R  of  fig. 4.  with  edges  F.v,  Tg,  I ’5  and  IV. 
This  region  R  represents  a  generic  rectangular  subdomain  in  the  domain  ft.  Let  nx  be  the  number 
of  gridpoints  in  T.y  (or  T*,-)  and  n2,  the  number  of  gridpoints  in  Tg  (or  ITv).  The  corner  points  are 
not  included  in  the  edges.  They  may  or  may  not  be  interior  to  ft. 
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For  t he  case  of  constant  coefficients  operators,  it  is  possible  to  describe,  in  terms  of  Fourier 
modes,  an  operator  Q  which  takes  boundary  values  on  one  of  the  edges  and  computes  the  solution 
on  the  gridpoints  adjacent  to  the  same  or  other  edge. 

Let  A  be  the  matrix  which  represents  the  discretization  of  the  differential  operator  in  ft.  If 
the  interface  r*.  where  k  =  .V,  S,  E  or  W,  is  interior  to  ft,  then  we  can  define  P the  submatrix 
of  A  that  represents  the  coupling  between  gridpoints  of  R  and  gridpoints  on  IV  Also  define  A/?, 
the  diagonal  block  corresponding  to  the  interior  points  of  R ,  in  other  words,  Ar  is  the  restriction 
of  the  differential  operator  to  the  region  R.  We  can  now  define  the  operator  Qki  which  represents 
the  interaction  between  the  edges  F*  and  T;  as: 


Qki  -  Pk  ar  Pi 


(5.1) 


For  constant  coefficients  operators,  when  T*  and  T/  are  parallel,  the  operator  Qki  is  diagonaliz- 
able  bv  Fourier  inodes.  For  example,  for  the  case  of  the  Poisson  equation  we  can  prove  the  rullowing 
lemma.  Here,  for  any  given  n,  \Vn  is  the  matrix  of  s  modes  of  dimension  n,  with  elements  given 
bv 


Wi 


■  sin 


IJTT 


(5.2) 


n  +  1  """  n  +  1 

for  i,j~  1, . . .,  n. 

Lemma  5.1.  Consider  the  Poisson  equation  on  a  domain  ft  which  contains  the  rectangular  region 
R.  Let  Qns  be  the  operator  that  represents  the  coupling  between  interfaces  IV  an d  Ts,  defined 
as  in  (5.1 ).  Then, 

WVQ/vsW'n,  =  Dns  , 

where  the  matrix  D.vs  is  diagonal,  with  diagonal  entries  given  by 


where 


Tj  = 


.ft) 


1  + 


n  +  m 


A 


and 


.(ft  - 


4  sin  ‘ 


2(tii  +  1) 

A  similar  expression  can  be  found  for  Qew- 
Also, 

Lftn,  QnnWtx,  =  DnN  1 

where  the  matrix  D,\y  is  diagonal,  with  diagonal  entries  given  by 


(5.3) 


(5.4) 


(5.5) 


djj  —  \/Yj 


1  ~  1? 


T!2  +  l 


(5.6) 


Similar  expressions  can  be  derived  forQss >  Qee  and  Qww- 

Proof.  Proofs  for  formulas  (5.3)  and  (5.6)  can  be  found  in  [3]  and  [4].  Here  we  give  a  different  - 
more  general  -  proof  using  direct  (or  tensor)  products.  The  matrices  Pt\  and  can  be  written  as: 


P.v  =  e[2)  0  /. 


(5.7) 
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r 


Ps  =  e‘2>  3  /, 


(5.8) 


where  //,  for  /  =  1.2,  is  the  identity  matrix  of  dimension  n/  and  is  the  3th  column  of  //. 
Tlie  matrix  .4/?  is  the  discrete  Laplacian  operator  on  the  region  R  and  it  has  the  following  block 
tridiagonal  form: 

/T  h  \ 

h  T 

,  h 

\  h  T } 


■4  R- 


(5.9) 


where  T  =  tridiag(  1,  —4, 1).  It  is  easy  to  prove  that 

VVniTWni  =  Dt  , 


where  Dj  =  diag(— 2  —  cr^ 1  ^ .  Then  we  have 


1 VniQNSWni  =  (e^&h) 


(  Dt  h 
h  Dt 

\ 


\ 

••  h 
h  dtJ 


-\ 


(5.10) 


By  reordering  the  equations  in  (5.10),  we  have: 

(Tl 


WniQNSWni  =  (I1Ze[2))T 


T2 


\ 

TnJ 


-X 


where  T}  =  tridiag(l,  -2  -  o]1',  1).  Therefore,  Wn^Q nsWni  is  diagonal  and  its  diagonal  elements 
are  given  by 

(2)T  j,_l  (2) 
el  en2  * 

which  can  be  proved  to  be  given  by  (5.3). 

Similarly,  we  can  prove  that  WniQ is  diagonal  and  its  diagonal  elements  are  given  by 

(2)rr-l  (2) 
cri2  j  »12  1 

which  can  be  proved  to  be  given  by  (5.6). 

I 

Operators  like  Qne,  on  the  other  hand,  which  represent  the  interaction  between  perpendicular 
edges,  are  not  diagonalizable  by  Fourier  modes.  Moreover,  they  are,  in  general  not  square,  but 
n\  by  ri2  rectangular  matrices.  It  is  possible,  however,  to  describe  the  elements  of  the  matrices 
Q.ve'^2  and  Wn,  <2,vivlT„2  for  constant  coefficients  cases. 

Lemma  5.2.  Consider  the  Poisson  equation  on  a  domain  Q  which  contains  the  rectangular  region 
R.  The  elements  of 

Qne  =  WrnQNEWn2 
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are  given  by 


q*E  = 


\/{n  i  +  1)(«2  +  1)  +  <7j 


„—ri  sin  -■‘•xt 

ni  -M _ n;  +  l 

(!)  ,  02) 


(5.11) 


for  i  -  1 . ti j  ami  j  =  1 . ;i2,  where  crj  =  4  sin2^^^ ,  for  l  =  1,2.  Similarly,  the  elements 

of  the  matrix 

Qnw  =  n2 

are  given  bv 


ti" 


2  sin  r^xr  sin 


nj+1 


^2  +  1 


y/{ni  +  1)(«2  +  U  <7,(1)  +  frj2) 

Proof.  The  eigenvalue  decomposition  of  the  matrix  Ar  is  well  known  and  it  is  given  by 

AR  =  (Wn2®Wni)  A  (U^  ©Un,) 
where  A  is  the  n\n-2  x  ni«2  diagonal  matrix  whose  diagonal  elements  are 

with  J  =  (j  -  l)ni  +  i.  lor  i  =  1, . . .,  «i  and  j  =  1, . . .,  712.  Also,  we  have 

P\v  =  h  0  • 


(5.12) 


(5.13) 


15.14) 


By  replacing  equations  (5.13)  to  (5.14)  in  (5.1)  and  then  applying  the  following  two  properties  of 
tensor  products: 

i)  (X  ®Y)T  =  XT  ®YTand 

ii)  (Xr  ©  y, )(X2  0  Y2)  =  (A'i-Y2)  0  Od^)  , 

we  have: 

Qnw  =  {(e\2)TW2)  0  Wx )  A"1  [W2  0  (5.15) 

and  therefore, 

Qnw  =  ((e< 2)TW2)  ®  h)  A"1  (/2  0  (W^))  (5.16) 

Then  we  can  see  that  the  j-th  column  of  (5.16)  is  given  by 


2  jit 

sin 


712  T  1  712  T 

from  which  (5.11)  follows. 

Similarly,  (5.12)  can  be  derived  by  using 


j  (Vj2)  /,  +  diag(<7,(1)))~1  W1  e(>) 


Pe  =  h  0  4  , 


(i) 


(5.17) 


instead  of  (5.7). 


6.  Appendix  B:  Proof  of  theorems  3.3  and  4.1 
Proof  of  Theorem  3.3 
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Proof. 

Theorem  3.3:  Let  the  aspect  ratio  lor  the  domain  ft2  in  fig-2  be  defined  as  a  =  ^^±1.  Then, 

a)  ||V'||i  <  A  0-733  and  ||V|U  <  0.733. 

b)  lll'T^l|2  <  !Iv'II:2  =  Ill~'ll2  <  ll^lli Halloo  <  0.54. 

c>  For  all  gridsizes  and  all  L-shaped  regions. 

AC(C4)  <  2.16  and  /C(C5)  <  2.16  .  (6.1) 

Proof,  (b)  follows  from  (a),  (c)  follows  from  (b)  and  from  equations  (3.24)  and  (3.25).  In  order  to 
prove  (a),  we  will  first  need  to  prove  some  lemmas  that  give  bounds  for  the  column  and  row  sums 
of  the  absolute  values  of  (3.28),  the  elements  of  the  matrix  V.  The  eigenvalues  (2.12)  of  A/4 
and  A/5  can  be  bounded  by 

l  /  \  v  .  •  3  ^  l  r.  o  \ 


Aj(n,mi,m2)  >  4sin 


2(n  +  l) 


for  all  j  =  1 .....  n  and 


\i(Tri2,n,nz)  >  4  sin 


2(m2  +  1) 


for  all  i  =  1, ... ,  m2.  It  is  easy  to  show  that 


1..  .  -  1  /(*i» 

K)l  £  vS  '■JTT 


where  the  function  /  is  defined  by 


f(x,y)  = 


'sinif  cosxf  ^/sinyf  cosj/f 
sin  2x~  +  sin  2y\ 


x,  =  ^2‘+  ;■  and  tjj  -  y^y.  Similarly,  we  have 

.  (6.6) 

2  m2  +  1 

The  column  and  row  sums  of  |v,y|  can  be  then  bounded  by  expressions  that  involve  the  integrals  of 
/  with  respect  to  x  or  with  respect  to  y.  The  following  lemma  gives  an  expression  for  the  integral 
of  /  with  respect  to  x,  for  a  fixed  y.  Since  f(x,y)  =  /(j/,x),  an  analogous  result  holds  for  the 
integral  of  /  with  respect  to  y,  for  a  fixed  x. 

Lemma  6.1.  Given  y  6  ( 0,1)  and  a,  b  €  (0, 1)  such  that  a  <  y  <  b, 
b 

/,  ,  2  7T  /  ,  .  ,  JT  .  7T .  ,.7T.7T\ 

/(x,  y)  dx  =  -j^-cosy-  -  g(smb-,sm  y-)  +  y(sma-,smy-)J  .  (6.7) 


where 


1  Z  +  y/2ziv  +  W  v/2  zw 

at  z.  w)  =  -  log - .  _ - +  arctan - 

J  2  h  z  -  +  w  -  -  w 
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Proof.  By  replacing  ;  =  sin  j:  ^  and  w  =  sin  in  (6.5)  and  defining 


F(z,w)  = 


\fzw 
z 2  +  w2 


we  get 


sin  b  - 


J  /(x,  y)  dx  =  ^  cos  y7^  J  F(z,w)dz 


z  4-  v2zw  -j-  w 


^COS4  [J2l0g  z  -  ^  F  W 
2 


+  arctan 


v/2l^\  !sin4? 


w  —  z 


=  ^Tcos^f  (x  _  5(sin6^,sin^^)  y(sina|,siny|)j 


(6.9) 


We  will  also  need  to  describe  t he  behiavor  of  f(x.y),  for  a  fixed  ^  £  (0,1),  in  the  interval 
x  £  (0, 1).  We  can  easily  see  that,  when  x,  y  £  (0, 1),  /(x,  y)  >  0.  In  the  next  lemma  we  prove  that 
/(.,;/)  has  one  and  only  one  relative  maximum  in  (0, 1). 

Lemma  6.2.  Given  y  £  (0, 1).  there  exists  a  unique  xm(y)  £  (n.  1)  such  that: 


ma x  f(x,y)  =  f(x',y)  , 

0<x<l 

/(.,  y)  is  monotonicaily  increasing  on  the  interval  (0,x”)  and  /(.,  y)  is  monotonically  decreasing  on 
(x",  1).  Moreover,  f{x',y )  is  bounded  by 


f(x\y)  <^=  cot  . 

Proof.  The  partial  derivative  of  /  with  respect  to  x  is  given  by: 

d  f  ,  .  2  tr  , .  .  2  tr 

—  =  ^x,j/)(sm  x-  -  r-Ksm^x-  -  r+)  , 

where  £(x,  y)  >  0  for  all  x,  y  £  (0, 1)  and 


(6.10) 


z±  =  |(1  +  sin2y|)±  ^(l  +  sin2t/| 


)2  -  sin2t/- 


It  can  be  shown  that 


sin  2y% 

0  <  z_  <  - —  <  1 


(6.11) 


and  z+  >  1.  Therefore,  §£  >  0  for  x  <  x'  and  <  0  for  x  >  x*,  where  x*  is  the  unique  solution 
in  (0. 1)  to 


sin  x  —  =  x_ 
2 


(6.11 


Therefore,  /  has  a  unique  maximum  in  (0, 1)  at  x*.  Moreover,  since  for  all  x  £  (0, 1), 


.  7T  .  7T 

/( x, .(/)  <  F(sin  x  — .  sin  y  — )  , 
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j 


where  F  is  defined  by  (6.9),  we  have 


r 


f(x',y)  <  max  F(z,  sin  y-)  =  ^  cot  y- 


We  can  now  prove  (a)  in  theorem  3.3: 

We  will  only  prove  the  inequality  [| V^jji  <  y/a  0.733.  The  proof  of  H^Hoo  <  0.733  is 

completely  analogous,  by  using  (6.4)  instead  of  (6.6). 

By  (6.6),  we  have 


7712 


2  \f0L  1 

=  max  >  Itml  <  - — — - -  max  >  f(xi,y;) 


i=l 


(C.  13) 


i=i 


Let  h  =  m--+1  =  xi .  Since  /  >  0  for  x,  y  6  (0, 1)  and,  by  lemma  6.2,  /  is  monotonic  in  the  intervals 
(0, x’(yj))  and  (xm(yj),  1),  it  is  easy  to  see,  by  using  graphical  arguments,  that 


m2  v 

f{x„yj)  <  f{x,Vj)dx  +  hf(x*(yj),y.:)  . 

i=i  { 


tfi  14) 


when  h  <  xM(yj).  On  the  other  hand,  when  h  >  x*(yj),  all  the  values  of  x;,t  =  l,...,m2  are  on 
the  interval  (x*.  1),  where  /  is  monotonically  decreasing.  Then,  we  have 


re?  m;  * 

hYL^Xi'yi>>  =  hf(h^yj)  +  h^f{xi,yi)  <  /  /(*, 2/j)dx  + 

1=1  i= 2  { 


hf(h,ijj )  .  (6.15) 


Let  us  first  assume  that  h  <  x‘{xjj).  By  (6.11)  and  (6.12),  we  can  prove  that  x"(y)  <  y  for  all 
y  E  (0, 1).  Then,  by  (6.7),  we  have: 


1 

J  f(x,yj)  dx  <  (tt  +  <7(sin/i|,sinyj|)) 


(6.16) 


because  y(z,  w)  >  0  for  w  <  z.  Define  the  function  G(h,(3)  =  #(sin /iy, sin/3fi j).  Then,  the  right 
hand  side  of  (6.16)  can  be  written  as: 


■(x  +  G(h,0))  , 


(6.17) 


with  3  =  ^  >  1.  By  differentiating  G  with  respect  to  3  we  can  see  that  '^4  <  0  for  all  h.  E  (0, 1) 


and  3  E  [1,  +oo).  Therefore,  G  decreases  with  /?,  i.e. 


G( h,  3 )  <  G ( h ,  1 )  -  lim  g(  z ,  w )  =  i  log  2  +  r  -  - 

*2+  2  2  —  v  2  2 


(6.18) 
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for  all  h  €  (0.  1)  and  J  E  [l.+oc).  Up  can  then  bound  (6.16)  by 


f  if  2 +  /2 


(6.19) 


On  the  other  hand,  by  (6.10).  we  have 

hf{x‘(’jj),yj)  <  -^=hcot 3li^  .  (6.20 

The  right  hand  side  of  (6.20)  can  also  be  proven  to  decrease  with  J  and  h  and  therefore  we  have 

hf(x*{>Jj\ >Jj)  <  j^=/icot/i|  <  —< ■  (6-21 

By  replacing  (6.19)  and  (6.21)  in  (6.14),  we  have 


=  1.4666 


and  therefore,  by  (6.13),  we  have 


!|T/||i  <  s/a  0.7333 


when  h  <  x'(yj). 

When  h  >  xn(y}),  by  (6.15)  we  have 


mi  j. 

J  /(x'J'i)dx+  hf{h,yj)  . 


(6.23) 


(6.24) 


By  (6.7), 


1 

J  f(x,yj)dx  <  cos  yj ^  (V  +  <7(sin.r'^,smt/;|)) 

•  (yj ) 


and.  by  (6.18),  we  have 


2  +  y/2 

2-/2 


(6.25) 


Since  }{h,  y.)  <  /(.r*( (/_,),  t/j ),  by  (6.21)  we  have 


hf(h,y/  < 


(6.26) 
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By  replacing  (6.25)  and  (6.26)  in  (6.24),  we  have 


m.2 

h^2f(x,,yj)  <  1.4666 

>'=  i 

and  therefore,  by  (6.13),  we  have 

II V Hi  <  x/aO.7333 

when  h  >  x By  (6.23)  and  (6.27).  we  proved  that  (6.23)  holds  for  all  h  <  1. 


(6-27) 


Proof  of  Theorem  4.1 

Proof.  Define  the  function 

We  can  easily  prove  that 


f{x)  = 


1  +  X  —  y/x 
1  —  2 


/(*)>  0.866  for  all  are  [0,1)  . 


By  (4.4)  and  (2.12),  we  have 
\j(n,  mA,  rn2)  -  (£_,(«,  m2)(  = 


„"»2  +  l 


2 


i  +  1  +  7j 7  j _ 

v^i  +  t  i  _  ^m-2  +  1  j  _  ^,m2  +  l 


1-7"11  1  -7* 


Since  -yj  <  1  and  m.\  <  m2,  we  have 


\j(n.  m,.  m2)  -  |fy(n,  m2)\  >  2 


1+7 


m2  +  l 


7 


m2  +  1 
2 


2  _  7m2  +  x  2  -  y 7,12  +  1 

.  '  J  '  ] 


=  2\/^  +  t/(7  r+1) 


>  1-73 \/ a,  +  -f 


I  cr ,  + 


CT* 


I  vj  +  ^r 


(6.28) 


(6.29) 


(6.30) 


Expressions  for  the  elements  of  Qse  and  Q,V£  of  (4.11)  are  given  in  appendix  A.  We  can  easily 
verify  that  the  elements  of  both  matrices  have  the  same  absolute  values.  Also,  both  A/8  and  A/g 
have  eigenvalues  that  are  bounded  from  below  by  4  sin  — & — 

J  "  2(m}  -H  )  ‘ 

By  (4.10),  (4.12)  and  (6.30),  we  can  see  that  ||V’||1  is  bounded  by 


ii^ih  < 


y/a  1 


\/0.866  l  2  m\  -f  1  i<j< 


-  max  V  V/sinj-I  cosar.f  ^/sint/jf  cos^-f 
1  i<j<n^  sin  2x,  |  +  sin  2yj\ 


(6.31) 


where  .r,  _  i/(m,  +  1)  and  yj  =  j/(n+  1),  for  i  =  l,...m,  and  j  =  l....n.  The  proof  of  theorem 
3.3  applies  now  to  the  expression  in  parenthesis. 
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